In geanders/improve_repro: Does not matter.

knitr::opts_chunk$set(echo = FALSE)

## Load packages with additional functionality beyond base R 
library(readxl)
library(tidyverse)
library(knitr)
library(purrr)
library(hms)
library(gridExtra)

# Read in data from the Excel spreadsheet template used to collect the
# data in the laboratory
growth_data <- read_excel("growth curve data_GR.xls",
                          sheet = "simplified_template")

# For each measurement time, calculate the time since the start of the
# experiment
growth_data <- growth_data %>% 
  mutate(sampling_delta_time = difftime(sampling_date_time, 
                                        first(sampling_date_time),
                                        units = "hours"))

# Convert to a format that will be easier to work with in plotting and
# data analysis (one row for each combination of measurement time and 
# test tube).
growth_data_tidy <- pivot_longer(growth_data, 
                                 -c(sampling_date_time, 
                                    sampling_delta_time),
                           names_to = "test_tube_id", 
                           values_to = "optical_density")

# From the test tube IDs, extract information like the growth condition (low 
# oxygen versus aerated) and a short tube ID (e.g., A1, L2). You can use
# regular expression functions (from the `stringr` package in the `tidyverse`)
# to do this.
growth_data_tidy <- growth_data_tidy %>% 
  mutate(test_tube_number = str_extract(test_tube_id, "[0-9]"),
         growth_condition = str_extract(test_tube_id, "[a-z_]*"), 
         growth_condition = str_replace(growth_condition, "_", " "), 
         growth_condition = str_to_title(growth_condition),
         short_tube_id = str_c(str_sub(growth_condition, 1, 1), test_tube_number)) %>% 
  select(sampling_date_time, sampling_delta_time, short_tube_id, growth_condition, optical_density)

ggplot(growth_data_tidy, 
       # Plot time since start on the x-axis, optical density on the y-axis,
       # and use color to show the growth condition (aerated versus low oxygen)
       aes(x = sampling_delta_time, 
           y = optical_density, 
           color = growth_condition)) + 
  # Add a line with different patterns for each test tube
  geom_line(aes(linetype = short_tube_id)) + 
  # Add points for each measurement
  geom_point() + 
  # Customize the labels for the legends and x- and y-axes
  labs(x = "Time (hours)",
       y = "Optical density at 600 nm",
       color = "Growth condition", 
       linetype = "Test tube number") + 
  # Customize the plot appearance
  theme_classic() + 
  # Add a title and subtitle
  ggtitle("Growth curve", subtitle = "Linear scale on y-axis") + 
  # Ensure that the y-axis ranges from 0 to 1
  ylim(c(0, 1))

ggplot(growth_data_tidy, 
       # Plot time since start on the x-axis, optical density on the y-axis,
       # and use color to show the growth condition (aerated versus low oxygen)
       aes(x = sampling_delta_time, 
           y = optical_density, 
           color = growth_condition)) + 
  # Add a line with different patterns for each test tube
  geom_line(aes(linetype = short_tube_id)) + 
  # Add points for each measurement
  geom_point() + 
  # Customize the labels for the legends and x- and y-axes
  labs(x = "Time (hours)",
       y = "Optical density at 600 nm",
       color = "Growth condition", 
       linetype = "Test tube number") + 
  # Customize the plot appearance
  theme_classic() + 
  # Add a title and subtitle
  ggtitle("Growth curve", subtitle = "Log-10 scale on y-axis") + 
  # Use a log scale for the y-axis. Ensure that the y-axis ranges 
  # from 0.01 to 1 (you can't start at 0, since that's undefined 
  # as a log transform)
  scale_y_log10(limits = c(0.01, 1))

\newpage

# Write a short function that can find the measurement time that is closest to, 
# but does not exceed, a specified time. For example, if you want to measure
# growth rate between 24 hours and 48 hours, you can use this function to find
# the measurement times closest to (while still being below) 24 and 48 hours.
find_closest_time <- function(data, check_time) {
  data %>%
    filter(sampling_delta_time < check_time) %>% 
    mutate(diff_from_check = check_time - sampling_delta_time) %>% 
    slice_min(n = 1, order_by = diff_from_check) %>% 
    pull(sampling_delta_time)
}

# Use the function to find the sampling times closest to (but not 
# over) 24 and 65 hours
closest_to_24 <- find_closest_time(growth_data, 24)
closest_to_65 <- find_closest_time(growth_data, 65)

# Calculate doubling times from the data
doubling_times <- growth_data_tidy %>% 
  filter(sampling_delta_time == closest_to_24 |
           sampling_delta_time == closest_to_65) %>% 
  group_by(short_tube_id, growth_condition) %>% 
  summarize(delta_time = diff(sampling_delta_time), 
            delta_optical_density = diff(log2(optical_density))) %>% 
  mutate(generation_time = delta_time / delta_optical_density)

a <- ggplot(growth_data_tidy, 
       aes(x = sampling_delta_time, 
           y = optical_density, 
           color = growth_condition)) + 
  geom_rect(xmin = closest_to_24, 
            xmax = closest_to_65, 
            ymin = -Inf, 
            ymax = Inf, 
            color = "beige", 
            fill = "beige") + 
  geom_line(aes(linetype = short_tube_id)) + 
  geom_point() + 
  labs(x = "Time (hours)",
       y = "Optical density at 600 nm)",
       color = "Condition type", 
       linetype = "Condition number") + 
  theme_classic() + 
  scale_y_log10() + 
  theme(legend.position = "none") + 
  ggtitle("Growth curves", subtitle = "Log scale") 

b <- ggplot(doubling_times,
            aes(x = short_tube_id, 
                y = generation_time, 
                fill = growth_condition)) + 
  geom_col() + 
  theme_classic() + 
  labs(y = "Doubling time (hours)", 
       x = "Test Tube ID", 
       fill = "Growth conditions") + 
  theme(legend.position = "top") + 
  ggtitle("Doubling times", subtitle = "24 to 65 hours")

grid.arrange(a, b, ncol = 2)

Table of estimated doubling times

# Create a table giving the estimated doubling times
doubling_times %>% 
  select(short_tube_id, growth_condition, generation_time) %>% 
  kable(col.names = c("Tube ID", "Growth conditions", "Estimated doubling time"), 
        digits = 2)